load libraries

library(tidyverse)
library(psych)
library(knitr)

set color palettes

palette = wesanderson::wes_palette("Zissou1", n = 4, type = "continuous")
palette_sub = wesanderson::wes_palette("Zissou1", n = 44, type = "continuous")

load data

# define variables and paths
sub_dir = "~/Dropbox (PfeiBer Lab)/FreshmanProject/Tasks/ROC-C/output/Himmel3_pilot/"
sub_pattern = "HL[0-9]{3}"
subjects = list.files(sub_dir, pattern = sub_pattern)
runs = c("run1", "run2", "run3")

# initialize data frame
data = data.frame()

# load data
for (sub in subjects) {
  for (run in runs) {
    file = paste0(sub_dir, sub, '/', sub, '_', run,'.csv')
    if (sub %in% c("HL001","HL002","HL003","HL004","HL005","HL006","HL007",
                   "HL008","HL009","HL010","HL011","HL012","HL013","HL014",
                   "HL015","HL016","HL017","HL018","HL019","HL020","HL021",
                   "HL022","HL023","HL024","HL025","HL026","HL027","HL028",
                   "HL029","HL030","HL031","HL032","HL033","HL034","HL035",
                   "HL036","HL037")){
          tmp = tryCatch(read.csv(file) %>% 
                     mutate(subjectID = sub,
                            run = run,
                            respCue = as.integer(as.character(respCue)),
                            respRating = as.integer(as.character(respRating)),
                            respEffort = as.integer(as.character(respEffort)),
                            batch = '1'), error = function(e) NULL)
          data = rbind(data, tmp)
    } else {
      tmp = tryCatch(read.csv(file) %>% 
                     mutate(subjectID = sub,
                            run = run,
                            respCue = as.integer(as.character(respCue)),
                            respRating = as.integer(as.character(respRating)),
                            respEffort = as.integer(as.character(respEffort)),
                            batch = '2') %>%
                       select(-craved), error = function(e) NULL)
      data = bind_rows(data, tmp)
      rm(tmp)
    }
  }
}

tidy data

data1 = data %>%
  #filter(subjectID %in% c("HL037","HL038","HL039","HL040","HL041","HL042","HL043")) %>%
  mutate(rtCue = ifelse(rtCue == "NaN", NA, rtCue),
         rtRating = ifelse(rtRating == "NaN", NA, rtRating),
         rtEffort = ifelse(rtEffort == "NaN", NA, rtEffort),
         action = ifelse(respCue == 1, "look",
                  ifelse(respCue == 2, "regulate", NA)),
         action = ifelse(cond == "LOOK" & is.na(respCue), "look",
                  ifelse(cond == "REGULATE" & is.na(respCue), "regulate", action)),
         choice = ifelse(cond %in% c("LOOK", "REGULATE"), "no",
                  ifelse(cond == "CHOOSE", "yes", NA)),
         ratingDiff = likingRating - respRating) %>%
  extract(col = foodPic, into = "category", regex = "(.*)[0-9]{2}", remove = FALSE) %>%
  group_by(subjectID) %>%
  mutate(trial = row_number()) %>%
  select(subjectID, run, action, choice, cond, respCue, trial, everything())

percent regulate

There are a few outliers, but most people chose to regulate around 40% of the time.

percent.plot = data1 %>%
  filter(choice == "yes") %>%
  group_by(subjectID, action) %>%
  summarize(n = n()) %>%
  spread(action, n) %>%
  mutate(percent.regulate = regulate / (look + regulate)) %>%
  arrange(percent.regulate)

ggplot(percent.plot, aes("", percent.regulate)) +
  geom_boxplot() +
  labs(x = "", y = "percent chosen to regulate")

percent.plot %>%
  ungroup() %>%
  select(percent.regulate) %>%
  as.data.frame() %>%
  describe() %>%
  select(n, mean, sd, min, max) %>%
  kable(digits = 2, format = "pandoc", caption = "group average")
group average
n mean sd min max
X1 44 0.47 0.09 0.25 0.67 
percent.plot %>%
  kable(digits = 2, format = "pandoc", caption = "by subject")
by subject
subjectID look regulate percent.regulate
HL006 36 12 6 0.25
HL007 18 6 12 0.25
HL003 24 12 18 0.33
HL044 36 18 NA 0.33
HL021 33 18 3 0.35
HL037 33 18 3 0.35
HL002 33 21 NA 0.39
HL011 33 21 NA 0.39
HL026 33 21 NA 0.39
HL042 33 21 NA 0.39
HL034 30 21 3 0.41
HL029 27 21 6 0.44
HL001 30 24 NA 0.44
HL005 30 24 NA 0.44
HL017 30 24 NA 0.44
HL018 30 24 NA 0.44
HL036 30 24 NA 0.44
HL041 30 24 NA 0.44
HL008 27 24 3 0.47
HL013 27 24 3 0.47
HL016 27 24 3 0.47
HL027 27 24 3 0.47
HL039 27 24 3 0.47
HL004 27 27 NA 0.50
HL009 27 27 NA 0.50
HL012 27 27 NA 0.50
HL014 27 27 NA 0.50
HL015 24 24 6 0.50
HL023 27 27 NA 0.50
HL025 27 27 NA 0.50
HL043 27 27 NA 0.50
HL020 24 27 3 0.53
HL022 24 27 3 0.53
HL010 24 30 NA 0.56
HL032 24 30 NA 0.56
HL033 24 30 NA 0.56
HL035 24 30 NA 0.56
HL038 24 30 NA 0.56
HL040 24 30 NA 0.56
HL019 21 27 6 0.56
HL028 21 27 6 0.56
HL031 21 27 6 0.56
HL024 18 33 3 0.65
HL030 18 36 NA 0.67

craving ratings

mean ratings by condition

looks like there is a main effect of action, doesn’t look like an interaction with choice

data1 %>%
  group_by(action, choice) %>%
  summarize(meanRating = mean(respRating, na.rm = TRUE),
            n = n()) %>%
  kable(digits = 2, format = "pandoc", caption = "mean craving ratings by condition")
mean craving ratings by condition
action choice meanRating n
look no 2.88 789
look yes 2.82 1188
regulate no 2.23 783
regulate yes 2.16 1071
NA yes 2.62 99 
plot.rating = data1 %>%
  group_by(subjectID, action, choice, batch) %>%
  summarize(meanRating = mean(respRating, na.rm = TRUE)) %>%
  filter(!is.na(action))

ggplot(plot.rating, aes(action, meanRating, shape = batch)) +
  geom_point(aes(group = subjectID), fill = "black", alpha = .05, size = 3) + 
  geom_line(aes(group = subjectID), color = "black", alpha = .05, size = 1) + 
  stat_summary(aes(color = action), fun.data = "mean_cl_boot", size = 1.5) +
  facet_grid(~choice) +
  scale_color_manual(values = c(palette[2],palette[4])) +
  labs(y = "mean craving rating")

ggplot(plot.rating, aes(choice, meanRating, shape = batch)) +
  geom_point(aes(group = subjectID), fill = "black", alpha = .05, size = 3) + 
  geom_line(aes(group = subjectID), color = "black", alpha = .05, size = 1) + 
  stat_summary(aes(color = choice), fun.data = "mean_cl_boot", size = 1.5) +
  facet_grid(~action) +
  scale_color_manual(values = c(palette[1],palette[3])) +
  labs(y = "mean craving rating")

individual differences

ggplot(plot.rating, aes(action, meanRating, color = choice, shape = batch)) +
  geom_point(aes(group =  interaction(subjectID, choice), color = choice)) +
  geom_line(aes(group = interaction(subjectID, choice), color = choice)) + 
  facet_wrap(~subjectID, ncol = 5) +
  scale_color_manual(values = c(palette[1],palette[3])) +
  labs(y = "mean craving rating")

difference in pre-task and task craving ratings by condition

data1 %>%
  group_by(action, choice) %>%
  summarize(meanRating = mean(ratingDiff, na.rm = TRUE),
            n = n()) %>%
  kable(digits = 2, format = "pandoc", caption = "mean difference in craving ratings by condition")
mean difference in craving ratings by condition
action choice meanRating n
look no 0.30 789
look yes 0.33 1188
regulate no 0.89 783
regulate yes 0.98 1071
NA yes 0.60 99 
plot.rating = data1 %>%
  group_by(subjectID, action, choice) %>%
  summarize(meanRating = mean(ratingDiff, na.rm = TRUE)) %>%
  filter(!is.na(action))

ggplot(plot.rating, aes(action, meanRating)) +
  geom_point(aes(group = subjectID), fill = "black", alpha = .05, size = 3) + 
  geom_line(aes(group = subjectID), color = "black", alpha = .15, size = 1) + 
  stat_summary(aes(color = action), fun.data = "mean_cl_boot", size = 1.5) +
  facet_grid(~choice) +
  scale_color_manual(values = c(palette[2],palette[4])) +
  labs(y = "difference in craving rating")

ggplot(plot.rating, aes(choice, meanRating)) +
  geom_point(aes(group = subjectID), fill = "black", alpha = .05, size = 3) + 
  geom_line(aes(group = subjectID), color = "black", alpha = .15, size = 1) + 
  stat_summary(aes(color = choice), fun.data = "mean_cl_boot", size = 1.5) +
  facet_grid(~action) +
  scale_color_manual(values = c(palette[1],palette[3])) +
  labs(y = "difference in craving rating")

craving rating means by subject and action (pre and task)

There seem to be habituation effects with pre-task ratings being higher than task ratings

plot.rating = data1 %>%
  group_by(subjectID, action) %>%
  mutate(sort = mean(respRating, na.rm = TRUE)) %>%
  filter(!is.na(action)) %>%
  gather(condition, rating, respRating, likingRating) %>%
  mutate(condition = ifelse(condition == "respRating", "task rating",
                     ifelse(condition == "likingRating", "pre rating", condition)))
  
ggplot(plot.rating, aes(reorder(subjectID, sort), rating, color = action)) +
  stat_summary(fun.y = mean, geom = "point", position = position_dodge(.9)) + 
  stat_summary(fun.data = mean_cl_normal, geom = "errorbar", position = position_dodge(.9), width = 0) +
  facet_grid(~condition) + 
  scale_color_manual(values = c(palette[2], palette[4])) +
  labs(x = "subject", y = "mean craving rating") + 
  theme(axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 80 rows containing non-finite values (stat_summary).

## Warning: Removed 80 rows containing non-finite values (stat_summary).

plot.rating = data1 %>%
  filter(!is.na(action)) %>%
  gather(condition, rating, respRating, likingRating) %>%
  mutate(condition = ifelse(condition == "respRating", "task rating",
                     ifelse(condition == "likingRating", "pre rating", condition))) %>%
  group_by(subjectID, action, choice, condition, batch) %>%
  summarize(meanRating = mean(rating, na.rm = TRUE))

ggplot(plot.rating, aes(action, meanRating, shape = batch)) +
  geom_point(aes(group = subjectID), color = "black", alpha = .05, size = 3) + 
  geom_line(aes(group = subjectID), color = "black", alpha = .05, size = 1) + 
  stat_summary(aes(color = action), fun.data = "mean_cl_boot", size = 1) +
  facet_grid(choice ~ condition) +
  scale_color_manual(values = c(palette[2],palette[4])) +
  labs(y = "mean craving rating")

mean craving ratings across trials

Looks like a slight decrease in craving ratings over time for both conditions, but it’s more pronounced in the regulate condition. Probably habituation, but we’ve also gotten feedback that the regulate mindset began to bleed over into look trials.

# linear trend
ggplot(filter(data1, !is.na(action)), aes(trial, respRating, color = action, linetype = choice)) + 
  geom_smooth(method='lm', se=FALSE, aes(group=interaction(subjectID, choice, action)), alpha=.06, size=.1) +
  geom_smooth(method='lm', se=FALSE, size = 1.5) + 
  scale_color_manual(values = c(palette[2],palette[4])) +
  labs(y = "mean craving rating")
## Warning: Removed 80 rows containing non-finite values (stat_smooth).

## Warning: Removed 80 rows containing non-finite values (stat_smooth).

# LOESS trend
ggplot(filter(data1, !is.na(action)), aes(trial, respRating, color = action, linetype = choice)) + 
  geom_smooth(method='loess', se=FALSE, aes(group=interaction(subjectID, choice, action)), alpha=.06, size=.1) +
  geom_smooth(method='loess', se=FALSE, size = 1.5) + 
  scale_color_manual(values = c(palette[2],palette[4])) +
  labs(y = "mean craving rating")
## Warning: Removed 80 rows containing non-finite values (stat_smooth).
## Warning in sqrt(sum.squares/one.delta): NaNs produced
## Warning: Removed 80 rows containing non-finite values (stat_smooth).

mean craving rating by run

Visualizing the decrease over time by collapsing across run

data1 %>%
  group_by(action, choice, run) %>%
  summarize(meanRating = mean(respRating, na.rm=TRUE),
            n = n()) %>%
  kable(digits = 2, format = "pandoc", caption = "mean craving ratings by condition and run")
mean craving ratings by condition and run
action choice run meanRating n
look no run1 2.89 264
look no run2 2.87 267
look no run3 2.87 258
look yes run1 2.82 390
look yes run2 2.84 393
look yes run3 2.80 405
regulate no run1 2.26 264
regulate no run2 2.23 261
regulate no run3 2.18 258
regulate yes run1 2.22 342
regulate yes run2 2.11 372
regulate yes run3 2.14 357
NA yes run1 2.66 60
NA yes run2 2.63 27
NA yes run3 2.42 12 
plot.rating = data1 %>%
  group_by(subjectID, action, choice, run, batch) %>%
  summarize(meanRating = mean(respRating, na.rm=TRUE)) %>%
  filter(!is.na(action))

ggplot(plot.rating, aes(run, meanRating)) +
  geom_point(aes(group=interaction(subjectID, action), color=action), alpha=.1, size=3) + 
  geom_line(aes(group=interaction(subjectID, action), color=action), alpha=.1, size=1) + 
  stat_summary(aes(color=action), fun.data = "mean_cl_boot", size = 1) +
  stat_summary(aes(group = action, color=action), fun.y = "mean", geom = "line", size = 1) +
  facet_grid(batch ~ choice) +
  scale_color_manual(values = c(palette[2],palette[4])) +
  labs(y = "mean craving rating")

liking ratings (pre-task)

ratings by category and frequency chosen

plot.rating = data1 %>%
  group_by(category) %>%
  summarize(meanRating = mean(likingRating, na.rm=TRUE),
            n = n(),
            sd = sd(likingRating, na.rm=TRUE),
            se = sd/sqrt(n),
            upper = meanRating + 2*se,
            lower = meanRating - 2*se)

ggplot(plot.rating, aes(reorder(category, -n), meanRating, color = category)) + 
  geom_point(aes(size = n)) + 
  geom_errorbar(aes(ymin = lower, ymax = upper), width = 0) + 
  # stat_summary(aes(group = category, color = category, size = n), fun.data = "mean_cl_boot", size = 1) +
  # stat_summary(aes(group = category, color = category), fun.y = "mean", geom = "line", size = 1) +
  scale_color_manual(values = wesanderson::wes_palette("Zissou1", 15, type = "continuous")) +
  labs(x = "category", y = "mean liking rating") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

distribution by subject and category

We might want to consider changing the image population script to only choose images >2 if possible. Right now it’s just choosing the bottom 20 images to maximize the spread between “highly desired” and “less desired” images, but since participants seem to be rating a number of images as 1s, we can change the script to prioritize more highly rated images first.

ggplot(data1, aes(likingRating, fill=category)) + 
  geom_histogram(stat = "bin", binwidth = .5) + 
  facet_wrap(~subjectID, ncol=5) +
  scale_fill_manual(values = wesanderson::wes_palette("Zissou1", 15, type = "continuous")) +
  labs(x = "liking rating")

effort ratings

effort rating means by subject and action

It’s more difficult to regulate than to look. Maybe a small interaction between action and choice.

data1 %>%
  group_by(action, choice) %>%
  summarize(meanEffort = mean(respEffort, na.rm=TRUE),
            n = n()) %>%
  kable(digits = 2, format = "pandoc", caption = "mean effort ratings by condition")
mean effort ratings by condition
action choice meanEffort n
look no 1.73 789
look yes 1.80 1188
regulate no 2.14 783
regulate yes 2.15 1071
NA yes 1.87 99 
plot.rating = data1 %>%
  group_by(subjectID, action, choice, batch) %>%
  summarize(meanEffort = mean(respEffort, na.rm=TRUE)) %>%
  filter(!is.na(action))

ggplot(plot.rating, aes(action, meanEffort, shape = batch)) +
  geom_point(aes(group = subjectID), fill = "black", alpha = .05, size = 3) + 
  geom_line(aes(group = subjectID), color = "black", alpha = .15, size = 1) + 
  stat_summary(aes(color = action), fun.data = "mean_cl_boot", size = 1.5) +
  facet_grid(~choice) +
  scale_color_manual(values = c(palette[1],palette[3])) +
  labs(x = "action", y = "effort rating")

ggplot(plot.rating, aes(choice, meanEffort, shape = batch)) +
  geom_point(aes(group = subjectID), fill = "black", alpha = .05, size = 3) + 
  geom_line(aes(group = subjectID), color = "black", alpha = .15, size = 1) + 
  stat_summary(aes(color = choice), fun.data = "mean_cl_boot", size = 1.5) +
  facet_grid(~action) +
  scale_color_manual(values = c(palette[1],palette[3])) +
  labs(x = "action", y = "effort rating")

individual differences

ggplot(plot.rating, aes(action, meanEffort, color = choice, shape = batch)) +
  geom_point(aes(group =  interaction(subjectID, choice), color = choice)) +
  geom_line(aes(group = interaction(subjectID, choice), color = choice)) + 
  facet_wrap(~subjectID, ncol = 5) +
  scale_color_manual(values = c(palette[1],palette[3])) +
  labs(y = "mean effort rating")

reaction times

Looks like some subjects might just be pressing buttons without really thinking (e.g. those around 500 ms). We might want to determine some cut off and throw out subjects under that time. On the whole, craving and effort ratings seem to be similarly difficult to answer (with a few subjects as notable exceptions).

RT boxplots by subject

Boxplots for craving rating and effort for each subject

plot.RT = data1 %>%
  gather("rating", "rt", rtEffort, rtRating) %>%
  group_by(subjectID, rating) %>%
  mutate(sort = median(rt, na.rm = TRUE))

ggplot(plot.RT, aes(reorder(subjectID, sort), rt, fill = rating)) +
  geom_boxplot() +
  scale_fill_manual(values = palette[1:2], labels = c("effort", "craving")) + 
  labs(x = "subject", y = "RT") +
  theme(axis.text.x = element_text(angle = 45,hjust = 1))
## Warning: Removed 221 rows containing non-finite values (stat_boxplot).

RT means by subject

Means and 95% CIs for craving rating and effort for each subject

plot.RT = data1 %>%
  gather("rating", "rt", rtEffort, rtRating) %>%
  group_by(subjectID, rating) %>%
  mutate(sort = mean(rt, na.rm = TRUE))

ggplot(plot.RT, aes(reorder(subjectID, sort), rt, color = rating)) +
  stat_summary(fun.y = mean, geom = "point", position = position_dodge(.9)) + 
  stat_summary(fun.data = mean_cl_normal, geom = "errorbar", position = position_dodge(.9), width = 0) +
  scale_color_manual(values = palette[1:2], labels = c("effort", "craving")) + 
  labs(x = "subject", y = "RT") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 221 rows containing non-finite values (stat_summary).

## Warning: Removed 221 rows containing non-finite values (stat_summary).